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One often-used approximation in the study of binary compact objects (i.e., black holes and neu- 
tron stars) in general relativity is the instantaneously circular orbit assumption. This approximation 
has been used extensively, from the calculation of innermost circular orbits to the construction of 
initial data for numerical relativity calculations. While this assumption is inconsistent with generic 
general relativistic astrophysical inspiral phenomena where the dissipative effects of gravitational ra- 
diation cause the separation of the compact objects to decrease in time, it is usually argued that the 
timescale of this dissipation is much longer than the orbital timescale so that the approximation of 
circular orbits is valid. Here, we quantitatively analyze this approximation using a post-Newtonian 
approach that includes terms up to order (Gm/(rc 2 )) 9 ^ 2 for non-spinning particles. By calculating 
the evolution of equal mass black hole / black hole binary systems starting with circular orbit config- 
urations and comparing them to the more astrophysically relevant quasicircular solutions, we show 
that a minimum initial separation corresponding to at least 6 (3.5) orbits before plunge is required 
in order to bound the detection event loss rate in gravitational wave detectors to < 5% (20%). In 
addition, we show that the detection event loss rate is > 95% for a range of initial separations that 
include all modern calculations of the innermost circular orbit (ICO). 

PACS numbers: 04.25.Dm, 04.30.Db, 04.40.Dg, 02.60.Cb 



I. INTRODUCTION 

The construction of astrophysically relevant initial 
data is a crucial prerequisite to using numerical relativ- 
ity as a predictive tool in gravitational wave astronomy, 
i.e. in producing gravitational waveform templates to 
be used in searches for signals in modern interferomet- 
ric gravitational wave detectors such as LIGO, VIRGO, 
GEO600, and TAMA300. After all, the comparison be- 
tween results from a numerical relativity evolution code 
and any observed astrophysical phenomena would be 
meaningless unless the initial data used in the numeri- 
cal evolution corresponds to an astrophysically realistic 
situation. In particular, numerical relativists are inter- 
ested in constructing astrophysically realistic initial data 
corresponding to the late stages of the inspiral process 
of compact objects (black holes and neutron stars), so 
that the transition from the final few orbits to the sub- 
sequent coalescence and formation of the final merged 
object can be numerically simulated using the full nonlin- 
ear theory of general relativity to accurately detail these 
highly nonlinear phenomena. To date, all initial data 
constructed for this purpose employ some sort of circu- 
lar orbit assumption in the construction of the initial 
data. For instance, many quasiequilibrium studies as- 
sume the existence of a helical Killing vector field in the 
construction of configurations that correspond to com- 
pact objiect binaries that are instantaneously in circular 
motion HH H 11 0,11 G3 IH 13 M Other 
methods employ turning point techniques to arrive at cir- 
cular orbit configurations 0, . It is typically argued 
that because the timescale of the gravitational radiation 
is longer than the orbital timescale of the binary, the as- 
sumption of quasiequilibrium is a good one. However, as 



the orbital separation of the compact objects decreases 
during the evolutionary progression of the binary, the 
timescale of the gravitational radiation process increases 
while the orbital timescale decreases. In other words, 
quasiequilibrium assumptions become increasingly inac- 
curate as the compact objects get closer. It is there- 
fore imperative to assess, in an astrophysically meaning- 
ful way, at what point during the inspiral process the 
quasiequilibrium approximation breaks down. In this pa- 
per, we present a method for calculating a lower bound 
on the error of the quasiequilibrium approximation for 
inspiraling binaries. This method is demonstrated by 
calculating this lower bound for equal mass black hole / 
black hole binaries and neutron star / neutron star bina- 
ries. 

One specific aspect of the quasiequilibrium approxima- 
tion is the circular orbit assumption, and it is precisely 
this aspect which we analyze in this paper. Our analy- 
sis is carried out within a post-Newtonian non-spinning 
point particle approximation. Using initial data corre- 
sponding to two particles in circular orbit, we solve the 
post-Newtonian equations of motion to find the resulting 
evolution of the binary system. The purpose of this cal- 
culation is twofold. First, we seek to quantitatively assess 
the circular orbit assumption in the construction of initial 
data in full general relativity. We do this by computing 
the correlation function of the gravitational wave signal 
produced from the evolution starting with particles in 
circular orbit with the gravitational wave signal from the 
true quasicircular evolution (i.e., the unique, "fully circu- 
larized" evolution of the binary system within the post- 
Newtonian approximation). This gives us a meaningful 
measure of the error introduced by the circular orbit as- 
sumption in the initial data. Second, we wish to have a 
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baseline with which to compare fully general relativistic 
calculations (i.e., calculations in numerical relativity) us- 
ing initially circular orbiting compact objects as initial 
data. 

We find that, as expected, the initial circular orbit as- 
sumption becomes progressively worse as the initial sepa- 
ration of the compact objects decreases. Specifically, we 
find that for equal mass black hole binaries, the initial 
separation must exceed a separation such that 6.3 (3.5) 
orbits remain before plunge in order to bound the de- 
tection event loss rate in gravitational wave detectors to 
< 5% (20%). This method is less efficient in calculating 
the lower bound of the error induced by the initial circu- 
lar orbit approximation for neutron stars, due to the fact 
that hydrodynamical effects will become important ear- 
lier in the evolution of the binary. As a result, the lower 
bound in the detection event loss rate for equal mass neu- 
tron stars that we calculate arc never larger than several 
percent. 

The remainder of the paper is organized as follows. In 
section m we describe the post-Newtonian equations of 
motion for two particles. We solve these equations for 
very large initial separation, thus producing the unique, 
"fully circularized" solution corresponding to the last 
1000 quasicircular orbits of the binary inspiral. In sec- 
tion IIIII wc formulate initial data corresponding to com- 
pact objects in initially circular orbits for various initial 
separations and compare the resulting gravitational wave 
signals to the true quasicircular gravitational wave sig- 
nals. In section IIVI we introduce a local definition of 
eccentricity which we find useful for comparisons to nu- 
merical relativity calculations, and compute the expected 
eccentricity as a function of initial separation that one 
can expect to find in numerical relativity calculations us- 
ing initial data corresponding to initially circular orbiting 
compact binaries. We compare these results with a fully 
consistent numerical relativity calculation of a neutron 
star binary inspiral using circular orbit initial data. In 
the conclusions, section we comment on the implica- 
tions our results have on the astrophysical relevance of 
current innermost circular orbit (ICO) calculations. 

II. QUASICIRCULAR BINARY EVOLUTION 

The general relativistic equations of motion for non- 
spinning point particles can be written in a post- 
Newtonian expansion as 
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where x — x 2 — x\ is the relative separation of the parti- 
cles, v — v 2 — Vi is the relative velocity between the parti- 
cles, r = \x\, n = x/r, m = mi+m 2 , r\ = m\m 2 /m 2 , and 



f = dr/dt. The post-Newtonian expansion is carried out 
in powers of e ~ m/r ~ v 2 (here, we have set G = c = 1). 

Of the non-radiative terms in the expansion (i.e., terms 
that are of integer order powers in the post-Newtonian 
expansion), only the 1PN and 2PN terms have been com- 
pletely determined 0, 0, 0] ; the 3PN terms have been 
calculated up to one numerical parameter [T^. l20l| . Of 
the radiative terms (e.g., n.5PN terms) in the expansion, 
only the 2.5PN El El [II and 3.5PN El terms have 
been completely determined. Employing an energy and 
angular momentum balance technique, the 4.5PN terms 
have been determined modulo 12 free "gauge" param- 
eters |2l|. Here, we numerically solve Eq. ^ using all 
post-Newtonian terms up to and including the first radi- 
ation reaction terms (i.e. up to and including the 2.5PN 
terms) for orbiting binaries. In order to obtain infor- 
mation regarding the error introduced by truncating the 
post-Newtonian expansion, we also calculate solutions 
which include the radiative 3.5PN and 4.5PN terms. We 
hereafter refer to this system of equations (i.e. up to and 
including the 2.5PN terms, plus the 3.5PN and 4.5PN 
terms) as the "4.5PN" equations, acknowledging that 
we are not including the conservative 3.0PN and 4.0PN 
terms. We have found empirically that the solutions to 
these 4.5PN equations of motion are highly insensitive to 
the 12 "gauge" parameters in the 4.5PN terms [2l]]. More 
specifically, the differences in the solutions to the 4.5PN 
equations of motion when the 12 parameters are ran- 
domly varied between —100 and 100 are orders of mag- 
nitude smaller than the differences in solutions obtained 
where the 4.5PN terms are dropped altogether. We take 
the binary to be in the x-y plane, which we coordinatize 
by the usual polar coordinates (r,4>)- Once the initial 
values of r, tf>, r, and (f> = dip/dt are set, the equations of 
motion completely specifies the solution. 

As shown in [22J, the effect of the radiation reaction 
terms in the post-Newtonian equations of motion for 
highly separated compact binaries (r > 1000 m) is to 
circularize the orbit of the binary. We refer to this state, 
and its subsequent solution to the post-Newtonian equa- 
tions of motion as the unique (up to rotations) "quasi- 
circular" solution. For an arbitrary binary separation 7*j, 
initial data which corresponds to this quasicircular orbit 
scenario is given at 2PN order (see |22j) as fi = with 
the square of the initial angular velocity <f>i satisfying the 
cubic equation 

= { r 2 tif-™{rU 2 f + (;Z-v)i-)\rUl)~ 

n n 
(15+^ + 27? 2 )(^) 5 . (2) 
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Since this initial data does not give a true quasicircular 
evolution when considering equations of motion of order 
2.5PN or higher, we use it as initial data for numeri- 
cal integrations starting with large initial separations T{ 
(typically r t > 100 m), and allow the equations of motion 
to fully circularize the orbit, after which we refer to the 
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FIG. 1: The orbital separation r as a function of the number 
of orbits for the quasicircular solution of the 2.5PN and 4.5PN 
equations of motion. We show the equal mass (mi = 7712) case 
along with the case of mass ratio mi/7712 = 10. The top panel 
displays the last 1000 orbits, the bottom panel displays the 
last 20 orbits. 



solution as the unique (up to spatial rotations) quasicir- 
cular solution. 

In Fig. 01 we plot the binary separation r as a function 
of cj> for the quasicircular solution to the 2.5PN and 4.5PN 
equations of motion for both the equal mass (mi — rn 2 ) 
case, as well as the case where the mass ratio mi/ma is 
10. We easily use enough resolution so that the trunca- 
tion error, if plotted as error bars on all curves in this 
paper, would be smaller than the thickness of the curves. 
We see that the 4.5PN solution displays a somewhat 
weaker radiation damping than the 2.5PN solution, in the 
sense that if one starts from a specific initial separation, 
there are more orbits until plunge for the 4.5PN solution 
as compared to the 2.5PN solution. For example, the bot- 
tom panel of Fig. ^ shows that for a circularized binary 
with orbital separation r = 15 m, the 2.5PN solution 
evolves for approximately 9 orbits until plunge (which, 
for definiteness, we define here as r = 2 m), whereas the 
4.5PN quasicircular solution evolves for approximately 
16 orbits until plunge. In Fig. |3 the ratio of the radial 
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FIG. 2: The ratio of the radial velocity |r| to the tangential 
velocity rcj> is plotted as a function of the binary separation r 
for the quasicircular solution for both the 2.5PN equation of 
motion and the 4.5PN equation of motion. We show the equal 
mass (mi = 7712) case along with the mass ratio mxjm-x = 10 
case. 

velocity \f\ to the tangential velocity rcf) as a function of 
separation r for the 2.5PN and 4.5PN quasicircular solu- 
tions is plotted, again for both the equal mass case and 
the mass ratio mi/ni2 = 10 case. In both cases, we see 
that for separations r > 4 m, the ratio of the radial ve- 
locity to the tangential velocity is smaller for the 4.5PN 
solution than for the 2.5PN solution; once the separation 
reaches 4 m, the merger is a fraction of an orbit away for 
both the 2.5PN and 4.5PN solutions. 



III. CIRCULAR INITIAL DATA AND 
SUBSEQUENT BINARY EVOLUTION 

The particular aspect of the quasiequilibrium approx- 
imation used in constructing initial data for binary sys- 
tems of compact objects in numerical relativity that we 
analyze here is the circular orbit assumption. In prac- 
tice, this assumption can take the form of an assumption 
of the existence of a helical Killing vector field; alterna- 
tively, turning point methods can be used to find approx- 
imate circular orbit configurations. The assumption of a 
circular orbit configuration as initial data in general rel- 
ativity, while fully consistent from a mathematical point 
of view, is certainly questionable from an astrophysical 
point of view. A casual inspection of Figs. and [21 re- 
veals that a circular orbit assumption in the initial data 
for a realistic binary system (which we effectively define 
as a fully circularized binary; we ignore here any capture 
scenarios and/or dynamically driven scenarios where the 
binary could have high eccentricities all the way down to 
the plunge phase) will be increasingly inaccurate as the 
initial separation decreases. We note that considerations 
of computational resources alone will limit accurate sim- 
ulations performed by numerical relativity codes in the 
coming decades to timescales of roughly 10 orbital peri- 
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FIG. 3: The angular velocity 4> = d<j>/dt is plotted as a func- 
tion of binary separation r for the quasicircular solution to 
the 2.5PN and 4.5PN equations of motion. Also shown is the 
angular velocity corresponding to instantaneously circular or- 
bits (f = and r = 0, see Eq. |1J at separation r. 
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FIG. 4: The evolved binary separation r as a function of time 
t is plotted for evolutions starting with circular initial data, 
with the initial separation fi = 12 m for the equal mass case 
(mi = m-2). Both 2.5PN solutions and 4.5PN solutions are 
plotted. For comparison, the quasicircular solution is also 
shown for both the 2.5PN and 4.5PN case. 



ods. From Fig. ^ this translates to initial binary sepa- 
rations ri < 20 m. The relevant question is therefore: at 
what initial orbital separation ri does the circular orbit 
approximation break down? In other words, will gravita- 
tional waveforms produced from an evolution code using 
initial data corresponding to compact objects in exactly 
circular orbit be, in some sense, close enough to wave- 
forms produced from the more astrophysically relevant 
quasicircular case? Here, we answer this question quan- 
titatively using the PN approximation described in the 
previous section. 

For an arbitrary initial separation r il we construct ini- 
tial data at time ti = corresponding to exact circular 
orbits by requiring the vanishing of the first and second 
time derivatives of the separation, fi = (dr/dt)\ t= o and 
fi = (d 2 r/dt 2 )\t=o- Assuming fi = 0, the expression for 
the instantaneous radial acceleration fi is given by the 
equations of motion, Eq. ^ as 

h = ntf + + A 1PN + A 2PN + A aPN + •••)• (3) 

We note that the radiation reaction terms (the n.5PN 
terms) do not enter in the prescription for circular 
data in the post-Newtonian approximation. Using the 
completely determined post-Newtonian Aipn and A 2 pn 
terms (see, e.g., ^H), the expression for the instanta- 
neous radial acceleration, Eq. can be written explicitly 
as 

h =nfa + %[ -1 -(l + 3»;)r?# + 2(2 + rj)^- 

r?(3 - 4rj)rt<J>t + ~f?(13 - Arfimr^ 2 - 
|(12 + 29t 7 )( — ) 2 . (4) 

We see that, to 2PN order, the equation r\ = (assum- 
ing fi is also zero) becomes a quadratic equation for the 



square of the initial angular velocity cf>i . The initial circu- 
lar orbit assumption, fi = and rt = 0, thus completely 
specifies the initial configuration for any arbitrary initial 
separation . The initial angular velocity </>j for this cir- 
cular condition is plotted in Fig.|3|as a function of initial 
separation Ti for the equal mass binary case. For compar- 
ison, the angular velocity for the 2.5PN and 4.5PN qua- 
sicircular solution is also shown. We see that the circular 
orbit approximation induces an artificial increase in the 
angular velocity as compared to the more astrophysically 
relevant quasicircular solutions. Differences between the 
two become apparent at roughly r = 6m, and by r = 3m 
the angular velocity for the circular orbit approximation 
is double that of the quasicircular solutions. 

In Fig. E| we plot solutions to the 2.5PN and 4.5PN 
equations of motion, using the initially circular condi- 
tions fi = and fi — as initial data, with the initial 
separation rj = 12m. Here, we show the equal mass case, 
mi = m,2. For comparison, the quasicircular solutions to 
the 2.5PN and 4.5PN equations of motion, where we have 
set t — when r = 12 m, are also shown. We see that 
the circular orbit approximation in the initial data in- 
duces oscillations (with a period of roughly 380 m) in the 
separation r, but that the evolutionary track of the two 
evolutions are similar. We see that for both the 2.5PN 
and 4.5PN circular initial data cases, the orbit is slightly 
eccentric, with t = corresponding to an apastron point 
in the orbit. Similar surveys of solutions starting with 
different initial binary separation parameter r, reveal the 
expected result that the oscillations in the evolved sepa- 
ration r using circular initial data decrease monotonically 
as the initial separation increases; i.e., the "error" in- 
duced by assuming circular initial data decreases with 
increasing initial binary separation rj. The fact that the 
circular orbit approximation induces an eccentricity in 
the orbit is not unexpected. It was shown in |23j | that 
allowing for eccentric orbits provided a better match be- 
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FIG. 5: The time averaged phase difference, A8 (Eq. 0, be- 
tween the quasicircular solution and the solution using in- 
stantaneously circular initial data for evolutions starting with 
initial binary separation n. The solid curve corresponds to 
solutions to the 2.5PN equations of motion and the dashed 
curve corresponds to solutions to the 4.5PN equations of mo- 
tion. 



tween post-Newtonian turning point calculation results 
and numerical relativity turning point calculation results. 

In order to get an idea of the phase errors induced by 
the circular orbit approximation, we calculate AO, de- 
fined by 



A6 = 



lj^dt(6 c (t)-6 gc (t)f 



f 



f 



(5) 



which is a measure of the time averaged phase difference 
between the two solutions to the post-Newtonian equa- 
tions of motion; 6 c {t) corresponding to the phase angle 
of the solution using circular orbit initial data and 9 qc (t) 
corresponding to the phase angle of the quasicircular so- 
lution, where each solution has an initial separation r% 
( AO is therefore an implicit function of initial separation 
ri). Here, we introduce a cutoff time tf, which we choose 
to be the time when either of the solutions reaches a bi- 
nary separation of r = 3m, since finite size effects become 
important for separations smaller than r — 3 m for black 
holes (recall that in harmonic coordinates, the horizon 
of a static non-rotating black hole is located at r = M, 
where M is the mass of the black hole). In Fig. [5] we 
plot A9 as a function of initial binary separation rj in 
the equal mass black hole case for both the 2.5PN and 
4.5PN solutions. We see that a maximum is reached for 
initial binary separations of r, = 8.8m and r, = 8.1m for 
the 2.5PN and 4.5PN solutions, respectively. These sep- 
arations correspond to roughly 1.5 orbital periods before 
the final plunge. 

From a gravitational wave detection point of view, the 
relevant measure of how "close" the solutions starting 
with circular initial data are to the more astrophysically 
realistic quasicircular solutions is determined by examin- 
ing the gravitational wave signal associated with each so- 
lution. We use the post-Newtonian formalism presented 
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FIG. 6: The "plus" (+) polarization gravitational waveform 
for coalescing, equal mass binaries. The top panel shows the 
radiation from solutions of the 2.5PN equations of motion and 
the bottom panel shows the radiation from solutions of the 
4.5PN equations of motion. Both panels show the waveform 
for the solution starting with circular initial data and the 
quasicircular solution. The initial separation in each case is 
r, = 10 m. The observation direction is in the orbital plane, 
with 6 = ty/2 and $ = 0. 



in |22j, |2J, |25j, |2|| , which give the two polarization states 
h+ and h x for our solutions as a function of observer dis- 
tance R and observation directions O and $. In Fig-EI we 
show gravitational waveforms for solutions to the post- 
Newtonian equations of motion using initial data starting 
at r, = 10 m. We see explicitly that the 2.5PN solu- 
tion plunges within approximately two orbits, while the 
4.5PN solution plunges within approximately four orbits 
(one orbital period is twice the period of the waveform) . 
We also see that the phase of the waveform correspond- 
ing to the circular initial data slightly lags the phase of 
the more astrophysically realistic quasicircular solution 
in both the 2.5PN and 4.5PN cases. The effect is less 
pronounced in the 4.5PN case. 

To quantify the difference between the gravitational 
waveform h c (t) obtained from circular initial data 
and the gravitational waveform h qc (t) obtained from 
the quasicircular solution, we define the correlation 
C[h c (t), h qc (t)] of the waveforms as 



C[h c (t),h qc (t)} 



Jo f h c (t)h qc (t - t) dt 



max ■ 



(h c (t)) 2 dt)(/ " (h qc (t-r)fdt) 



(6) 



Notice that the quantity inside the braces in the defini- 
tion of the correlation C[h c (t), h qc (t)] depends on both 
the lag time r and the cutoff time tf for the time inte- 
grations. In our case, we specify the cutoff time tf as the 
earliest time when either of the waveforms corresponds 
to some "final" binary separation rf. For black hole bi- 
naries, we expect tidal effects to become important as 
r < 3 m, so we set rf = 3 m. For neutron star equal 
mass binaries, hydrodynamical effects will become im- 
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FIG. 7: The event capture rate, C[h c (t), h qc (t)] (where 
C[h c (t),h qc (t)] is defined in Eq. |SJ, of 4.5PN gravitational 
waveforms using circular initial data as a function of initial 
binary separation r». The gravitational waveform observer is 
located in the orbital plane O = 7r/2, with observer polar 
angles $ = 0, tt/8, 2tt/8, 3tt/8, 4tt/8, 5tt/8, 6tt/8, and 7tt/8 
shown as thin dotted curves. The thick solid curve labeled 
"4.5PN" is the minimum correlation for all observation polar 
angles $. This same minimum is plotted for the 2.5PN case 
and is labeled "2.5PN". 



portant (depending on details of the equation of state) 
at around r = 9 m, in which case we set rj = 9 m. Once 
the cutoff time tf is determined, the quantity inside the 
braces in Eq. is a function of the lag time r; the fi- 
nal correlation C[h c (t), h qc (t)] is determined by varying 
the lag time r and finding the maximum value of the 
expression in braces. It is important to note that by in- 
troducing a cutoff of the gravitational waveform (which 
is necessary in a practical sense, due to the fact that the 
post-Newtonian approximation we are using breaks down 
when tidal effects become important), we are computing 
a lower bound of the error induced by using circular initial 
data for inspiral simulations; a more realistic calculation 
of the error would include all of the nonlinear phenomena 
occurring after tidal effects become important, which we 
are necessarily neglecting. 

In Fig. [7| we plot the cube of the correlation function 
C[h c (t),hq C (t)], Eq. El between the waveform computed 
using circular initial data and the waveform of the quasi- 
circular solution as a function of initial binary separation 
r,;, where we have set rj = 3 m. This gives us an esti- 
mate of the lower bound of the event loss rate induced by 
using circular initial data instead of the more astrophys- 
ically correct quasicircular initial data in the case where 
the compact objects are black holes. For each initial sep- 
aration r"i, we can find the minimum of the correlation 
over all observation polar angles <&, and take this to be 
the measure of the lower bound of the error induced by 
the circular orbit approximation in the initial data. This 
curve is labeled "4.5PN" in Fig. [7\ The oscillations in 
this curve, which have a period of slightly more than one 
additional orbit in the evolution of the binary, correspond 
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TABLE I: The minimum initial binary separation Vi for spe- 
cific event loss rates due to errors in the initial data introduced 
by the circular orbit approximation for black holes (77 = 3m), 
as predicted by the 2.5PN and 4.5PN equations of motion. 
In brackets are the number of orbits remaining until plunge, 
starting at this initial binary separation (see Figs. and 01. 



precisely to additional oscillations in the r vs. t curves 
(see Fig. |2J. That is, as the initial binary separation n 
increases, more oscillations in the r vs. t curve are al- 
lowed before the binary plunges, and the oscillations in 
Fig. [7| occur as a result. These oscillations are also evi- 
dent in the time averaged phase difference, A9 (Fig. EJ, 
and occur for the same reason. 

We see from Fig. [7] that, for any given initial separa- 
tion rj, the correlation between the waveform obtained 
from circular initial data and the waveform obtained from 
the quasicircular solution is larger for the 4.5PN solutions 
than for the 2.5PN solutions. However, recall from Fig.Q] 
that, for a given initial separation n, there is quite a large 
difference in the number of orbits before final plunge be- 
tween the 2.5PN case and the 4.5PN case. From the 
standpoint of numerical relativity, the number of orbits 
before plunge is perhaps the more important criterion in 
parameterizing the orbital separation of the binary. In 
Fig. |21 we recast the results in Fig.[7| plotting the cube of 
the correlation function C[h c (t),h qc (t)], Eq. El between 
the waveform computed using circular initial data and 
the waveform of the quasicircular solution as a function 
of the number of orbits until final plunge. We see that 
the profile of this event capture rate is similar for both 
the 2.5PN and 4.5PN cases, the major difference between 
the two being the phase of the oscillations. This gives us 
confidence that the post-Newtonian approximation has 
converged sufficiently to provide us with an accurate es- 
timate of the lower bound of the error induced by the 
circular orbit assumption in the construction of initial 
data. The results of Figs. [7| and 0D are summarized in 
tabular form in Table |U We can see that using circular 
initial data with configurations that have 2 orbits until 
plunge or less will result in a 50% or more event loss rate, 
and that one must start with configurations of n > 10 m 
(with more than 3.5 orbits until plunge) in order to have 
a chance of bounding the event loss rate to below 20%. 

For neutron stars, we expect tidal effects to become 
important at binary separations r of between 7 m and 
9 m, depending on the equation of state of the nuclear 
matter. As a result, the point particle approximation 
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FIG. 8: The event capture rate, C[h c (t), h qc (t)\ 3 (where 
C[hc(t),hqc{t)\ is defined in Eq.^J, of gravitational waveforms 
from inspiraling black holes (77 = 3 m) using circular initial 
data as a function of the number of orbits until plunge (i.e., 
until r = 2 m). The observer is located in the orbital plane 
O = 7r/2, and the event capture rate is minimized over all po- 
lar angles $ (see, e.g., Fig.0. Results using both the 4.5PN 
and 2.5PN equations of motion are shown. 
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FIG. 9: The event capture rate, C[h c (t), h qc {t)] 3 (where 
C[h c (t), h qc (t)] is defined in Eq.[5J, of gravitational waveforms 
from inspiraling neutron stars (77 =9 m) using circular ini- 
tial data as a function of the number of orbits until separation 
r = 77 —9 m, at which point hydrodynamics! effects will be- 
come important. The observer is located in the orbital plane 
<d = it/2, and the event capture rate is minimized over all po- 
lar angles $ (see, e.g., Fig. Q. Results using both the 4.5PN 
and 2.5PN equations of motion are shown. 



breaks down much earlier for neutron star models than 
for black hole models. This will necessarily result in a 
corresponding decrease in sensitivity when determining 
the errors induced by using circular orbit initial data. 
However, this does not necessarily imply that the errors 
induced by assuming initially circular orbits is less for 
neutron star binaries than those for black hole binaries 
at similar separations, but instead is simply an artifact 
that this particular model (spinlcss point particle model) 
is worse for neutron star binaries than for black hole bi- 
naries; after all, we are only computing lower bounds to 
these errors. 

In Fig. EH we set the cutoff separation 77 to 9 m (i.e., 
assuming neutron star binaries) and plot the cube of the 
correlation function C , Eq.|Sl between the waveform com- 
puted using circular initial data and the waveform of the 
quasicircular solution as a function of initial binary sep- 
aration r, ; . As expected, we observe a notable decrease 
in sensitivity to errors induced by assuming initially cir- 
cular data as compared to the black hole case (77 = 3 m, 
see Figs. 0and In the case of equal mass binary neu- 
tron stars, the lower bound on the event loss rate that 
we compute is lower than 2%. 



IV. CIRCULAR INITIAL DATA: COMPARING 
POST-NEWTONIAN AND FULL NUMERICAL 
RELATIVITY CALCULATIONS 

Numerical relativity has progressed to the point where 
stable, multiple orbit numerical simulations of compact 
objects are now possible. In |27|. we present fully general 
relativistic simulations of binary neutron stars using con- 
formally flat, quasiequilibrium initial data corresponding 



to an equal mass, corotating binary. The construction 
of the initial data, in addition to being a solution of the 
constraint equations of general relativity, assumed the ex- 
istence of a timelike helical Killing vector field. That is, 
the binary is assumed to be instantaneously in a circular 
orbit. We wish to compare the numerical relativity calcu- 
lations obtained in |27| with the post-Newtonian formu- 
lation from the previous sections. One difficulty in per- 
forming any such comparison is due to the fact that each 
calculation was performed in different coordinate sys- 
tems. The post-Newtonian calculations are performed in 
a harmonic gauge {g a ^T^ a (j = 0), whereas numerical rel- 
ativity simulations typically use a variety of other (both 
local and elliptic-type) gauge conditions to specify the co- 
ordinates. In particular, the relativistic calculations of bi- 
nary neutron stars presented in |27^ used a quasi-isotropic 
coordinate system; the initial spatial slice is conformally 
flat, but the subsequent evolution of the full Einstein 
field equations drives the 3-metric away from conformal 
flatness. Nevertheless, we can attempt a comparison be- 
tween the harmonic coordinate post-Newtonian calcula- 
tions with the quasi-isotropic coordinate numerical rel- 
ativity calculation by using a simple coordinate trans- 
formation that takes the single, stationary, non-rotating 
black hole solution of the Einstein equations written in 
harmonic coordinates |28j| to the same solution written in 
isotropic coordinates. Let (thar , fhar) and (ti SO ,ri SO ) be 
the (time, radial) coordinate pair for the harmonic and 
isotropic coordinatizations, respectively. They are re- 



NS-A (numerical relativity NS simulation) 
2.5PN, circular initial data 
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FIG. 10: The isotropic separation of binary evolution calcula- 
tions as a function of the number of orbits. Each calculation 
uses a circular orbit assumption as initial data. The solid 
curve refers to a fully consistent numerical relativity calcula- 
tion of binary neutron stars using conformally flat , quasiequi- 
librium initial data (NS-A from reference [27111 . The dotted 
and dashed curves correspond to 2.5PN and 4.5PN calcula- 
tions, respectively, starting with the same initial isotropic sep- 
aration as that of NS-A, assuming initially circular orbits. 
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(7) 



We use this relationship to map the harmonic coordi- 
nate separation of the post-Newtonian calculations done 
in harmonic coordinates to the isotropic separation of 
the fully general relativistic calculations in |27j done in 
(quasi-) isotropic coordinates. 

In Fig. 1101 we plot the isotropic, center of mass sepa- 
ration from the numerical relativity calculation NS-A of 
reference |27| as a function of the number of orbits. For 
comparison, we calculate the same quantity using the 
2.5PN and 4.5PN equations of motion from the previous 
sections, starting with the same isotropic binary separa- 
tion as that of NS-A and using the initially circular or- 
bit assumption. The frequency of the oscillations in the 
separation is slightly different for the numerical relativ- 
ity calculation NS-A as compared to the post-Newtonian 
calculations. There are many possible explanations for 
this. First of all, the calculations are done in different 
coordinate systems. While an attempt has been made to 
match the coordinates used to measure the binary sepa- 
ration (i.e., Eq.0| in the two cases, the differences in the 
angular and/or time coordinates could explain the differ- 
ence. Second, the numerical relativity calculation NS-A 
contains finite size effects, whereas the post-Newtonian 
approximation does not. Thirdly, the neutron stars in the 
initial data used for the NS-A numerical relativity calcu- 
lation contain an unphysically high amount of spin, since 
the co-rotation assumption was used in its construction. 
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FIG. 11: An example of the definition of An, 2, used in the 
definition of the generalized eccentricity, Eq. [5] Here, r(t) 
represents the 2.5PN solution using circular initial data at an 
initial binary separation of fj = 14 m. The general apastron 
points are denoted by circles, and ro(t) is the line connecting 
the generalized apastron points. 



The post-Newtonian approximation we are using does 
not take spin into account. Thus, any spin-orbit cou- 
pling in the NS-A calculation will not be present in the 
post-Newtonian approximation. Finally, the numerical 
errors in the calculation NS-A still could be large enough 
to explain the difference. 

However, the magnitudes of the oscillations in the sep- 
arations shown in Fig. 1101 which are a measure of the 
eccentricity of the orbits, are comparable. Notice that 
the separation r is a monotonically decreasing function 
of time for all three calculations presented in Fig. EH 
We note that this is universally true for post-Newtonian 
calculations using circular initial data. We therefore gen- 
eralize the concept of an apastron point (typically defined 
as a local maximum of separation r) for monotonically 
decreasing binary separation functions r(t) as follows: if 
r(t) is monotonically decreasing, then we define a gen- 
eralized apastron point to be those points where dr/dt 
attain a local maximum. We note that the numerical rel- 
ativity calculation NS-A, as well as the 2.5PN and 4.5PN 
solutions using circular initial data, are initially (t = 
corresponds to orbits in Fig. llOfl at generalized apas- 
tron points. Given any two consecutive generalized apas- 
tron points of a separation function r(t) at, say, times t\ 
and t2, we can define the eccentricity as follows. Let 
(r) 12 be the time averaged value of the binary separation 
r between the two generalized apastron points: 



tfdtrjt) 

(<2-tl) ' 



(8) 



Let r (t) be a linear function of time that intersects the 
points r(ii) and r(t<i). We define Ari,2 to be the range 
of r(t) — ro(t) within the interval t = [£1,^2] (see Fig. II IB . 



The generalized eccentricity e is then defined as 



The advantages of this particular definition of the gener- 
alized eccentricity are twofold. In the first place, Eq. [5] 
reduces to the usual definition of eccentricity in the New- 
tonian limit. Secondly, it is defined in terms of quantities 
local to the binary, as opposed to quantities defined at 
spatial infinity. This makes it particularly useful in nu- 
merical relativity, where invariant measures of the binary 
separation can be straightforwardly computed. 

In Fig. El we plot the generalized eccentricity e (Eq.[3J) 
for initially circular orbiting equal mass binaries as a 
function of initial binary separation n. The first gen- 
eralized apastron point used in each calculation is at the 
initial time t = 0. Shown are calculations using both 
the 2.5PN and 4.5PN equations of motion. At small val- 
ues of initial separation r-j (less than j-j = 10.07 m for the 
2.5PN equations of motion and rt = 9.18m for the 4.5PN 
equations of motion) the solutions r(t) have monotoni- 
cally decreasing first derivatives dr/dt, and thus do not 
exhibit multiple generalized apastron points (one could 
define the generalized apastron points in this case to be 
points where d 2 r/dt 2 obtain local maxima, thereby fur- 
ther generalizing the definition of apastron points, but 
we have not done so here). For comparison, we show 
as a solid triangle in Fig. ^| the results from numerical 
relativity simulation NS-A from reference [27| of equal 
mass binary neutron stars using conformally flat, quasi- 
equilibrium initial data. In order to make the compari- 
son more meaningful, we use Eq. [7| to relate the quasi- 
isotropic coordinate separation from the numerical rel- 
ativity calculation to the harmonic coordinates used by 
the post-Newtonian calculation (although the difference 
between the two are quite small for r.i — 11.6 m). As 
can be seen, the eccentricity observed in the numerical 
relativity calculation NS-A lies between the 2.5PN and 
4.5PN eccentricity predictions. Again, this comparison 
must be viewed as preliminary, as there still may exist 
significant numerical errors (i.e., boundary and/or trun- 
cation errors) in the numerical relativity simulations. 

V. CONCLUSIONS 
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FIG. 12: The generalized eccentricity e, Eq.0 as a function of 
initial separation n for initially circular orbiting equal mass 
binaries. Shown are results for both the 2.5PN and 4.5PN 
equations of motion. Also shown as the solid triangle is the 
eccentricity e of the numerical relativity simulation of binary 
neutron stars NS-A from reference |27jl . 



be required to perform at least 3.5 orbits (or more, de- 
pending on the accuracy required, see Table HJ) in order 
to produce astrophysically meaningful waveforms to be 
used as search templates in modern interferometric grav- 
itational wave detectors. 

Finally, we comment on the astrophysical relevance of 
recent innermost circular orbit (ICO) calculations (for a 
review, see OS HI). We note that all current calcula- 
tions of ICO configurations for equal mass black holes 
have an angular velocity m<j> of as low as 0.06 and as 
high as 0.18. From Fig. |3 we see that this range corre- 
sponds to binary separations of r = 6 m or less. Also, 
all current ICO calculations are based on the circular or- 
bit approximation. However, we have shown in this paper 
that the circular orbit approximation has completely bro- 
ken down for equal mass black holes with separations of 
r < 8.5 m. This calls into question the astrophysical sig- 
nificance of ICO calculations based on circular orbit as- 
sumptions. Certainly, any gravitational waveforms pro- 
duced from numerical relativity codes using initial data 
that corresponds to ICO configurations computed from 
circular orbit assumptions will be of little or no astro- 
physical significance. 



In this paper, we have performed an analysis of the 
circular orbit approximation which is often used in the 
analysis of binary compact objects in general relativ- 
ity. Using astrophysically relevant measures of the er- 
rors induced by the circular orbit approximation, we have 
shown that for equal mass binary black holes, the ap- 
proximation has completely broken down for binary sep- 
arations of r < 8.5 to, which correspond to separations 
where roughly 2 orbits or less remain before the final 
plunge. We conclude that numerical relativity evolution 
codes using circular orbit initial data configurations will 
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